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Abstract. In this paper, we introduce an ASEP-like transport model for bidirectional 
' motion of particles on a multi-lane lattice. The model is motivated by in vivo 

■ experiments on organelle motility along a microtubule (MT), consisting of thirteen 
protofilaments, where particles are propelled by molecular motors (dynein and kinesin). 
In the model, organelles (particles) can switch directions of motion due to "tug-of-war" 
events between counteracting motors. Collisions of particles on the same lane can be 
cleared by switching to adjacent protofilaments (lane changes). 

■ We analyze transport properties of the model with no-flux boundary conditions 
, at one end of a MT ("plus-end" or tip). We show that the ability of lane changes 

can affect the transport efficiency and the particle-direction change rate obtained 

■ from experiments is close to optimal in order to achieve efficient motor and organelle 
iy) \ transport in a living cell. In particular, we find a nonlinear scaling of the mean tip size 



(the number of particles accumulated at the tip) with injection rate and an associated 
phase transition leading to pulsing states characterized by periodic filling and emptying 
of the system. 



Keywords: traffic and crowd dynamics, molecular motors (theory), stochastic particle 
dynamics (theory), phase transition 

1. Introduction 

The cytoplasm of living cells contains a complex network of fibres (the cytoskeleton) 
that helps to maintain cell shape by providing structural support, and that facilitates 
intracellular transport. This transport is mediated by specialized mechano-proteins, 
the so-called molecular motors, that utilize ATP to move vesicles, organelles (or other 
cargo) in a particular direction along the cytoskeleton pQ. Microtubules (MTs) are one 
type of cytoskeletal element that consist of thirteen oriented protofilaments [2] (that 
we refer as lanes), each of which can support bidirectional transport. Long-distance 
transport along MTs is powered by kinesin and dynein, where kinesin takes its cargo to 
the polymerization-active plus-end of a MT while dynein walks towards the minus-end 
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Various lattice models have been proposed to understand the role of cooperation 
and competition between motors involved in such bi-directional transport. These 
are typically generalizations of asymmetric simple exclusion processes (ASEP) in one 
dimension [3j H] where motors are represented by particles that move on a lattice. 
Already for unidirectional motor transport models [HI El [7] there are nontrivial effects 
such as traffic jams due to the mutual exclusion of particles on sites of the MT [5j [8]. 

Bidirectional transport has previously been modelled by an exclusion process with 
binding and unbinding events [91 QUI El EE] or a lattice model with site sharing 
[13] . Change in direction of individual cargos may be due to "tug-of-war" [j] events 
[TH [151 [161 ITT] and this can be incorporated into ASEP models by assuming that motion 
of a single type of particle in different directions takes place on different lanes (and 
change of direction is modelled simply by change of lane) [HI HH [20] or that particles 
may change type [21]. In the latter case, counter- moving particles approaching on the 
same lane can only pass if at least one of them can change lane. The lane-change 
rules in the latter model are based on experimental observations that (a) dynein can 
switch between lanes on the MT at a certain rate [22], whereas (b) kinesin remains on 
a single lane [23], and (c) no long-lasting traffic jam is observed between particles away 
from the tip. However, recent reports show that kinesin can change lanes to overcome 
obstacles on the MT [23]. Thus the situation in the living cell is less clear and one aim 
of this paper is to investigate the impact of different lane-change rules on bidirectional 
transport. 

We study a generalization of the models from [20] I2l] for bidirectional motion 
motivated by in vivo experimental observations of bidirectional motion of dynein 
particles near the hyphal tip of Ustilago maydis. The particles represent dynein that 
is transported towards the plus-end of a MT by kinesin-1 and towards the minus-end 
under its own power. As the dynamics of the MT in this system are comparatively stable 
[25] and binding/unbinding events are rare in this system [26], we focus on modelling 
bidirectional transport of particles that remain bound at all times to a fixed section of 
MT where one end is a plus-end (with no-flux boundary conditions are applied) while 
at the other end we inject particles at a given injection rate. In particular we explore 
the behaviour of the tip accumulation (or tip size) as a function of this injection rate. 

The paper is organized as follows: In Section [2] we introduce a general multi-lane 
model for bidirectional motion and demonstrate that it includes, as special cases, the 
models of [201 12T] . Section [3] examines the influence of collisions and lane changes on 
bidirectional transport in this model. In particular we find cross-lane diffusion due to 
lane changes, and for low injection rates we find an approximately linear scaling of the tip 
size with injection rate as discussed in [20]. However, for higher injection rates we find 
there is a nonlinear growth in the tip size (depending on the lane-change rules) due to 
trapping of particles near the plus-end, and a singularity at a finite critical injection rate. 
We suggest that is associated with a phase transition of the system. We also study how 

J "tug-of-war" events refer to simultaneous and competitive activity of counteracting motors on the 
same organelle. 
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the lane change rules influence the tip structure and discuss measures of efficiency for 
the bidirectional transport. Interestingly, we find that the in vivo parameters obtained 
from living cell images [2T] are close to optimal for effective bidirectional transport, in 
that they balance efficient transport to the tip against the average delay at the tip. In 
Section H] we describe novel pulsing states (with whole-system approximately periodic 
filling and emptying) and filled states that can both appear at injection rates beyond 
the critical value. We finish with a discussion of biological relevance, implications and 
generalizations of these results in Section [51 

2. A multi-lane model for bidirectional transport 

Consider a bidirectional transport model where particles move on a lattice of spatial 
locations consisting of M adjacent lanes around the circumference of a cylinder J§| Each 
oriented lane is discretized into N sites between a plus-end and minus-end, as illustrated 
in Figure [TJ Particles are of two types; plus-type particles move towards the plus-end, 
while minus-type particles move towards the minus-end; a single particle represents a 
bound pair of opposite-directed motor proteins that is pulled along the lane by one of 
the motors being bound to the MT and the type of the particle corresponds to which 
of the motors is currently bound to the MT. 

We identify each location along the cylinder by a pair (/, i) where I G {1, • • • , M} 
denotes the lane and i G {1,---,N} the site along the lane and let r l ±i = 1 or 
represent the presence or absence of a plus-type (minus-type) particle at location 
Each location is occupied by at most one particle (i.e. t\ i + r l _ i can only be or 1) and 
particles move from one location to another at given rates. We also allow the possibility 
that each particle can change from one type to the other at rates representing the 
resolution of brief "tug-of-war" events [El [15j [16], [T7] between opposite oriented motor 
proteins bound to the particle. This change of type may or may not be associated with 
a change of lane. 

The plus- and minus- type particles may have transition rates that describe motion 
and particle-type change and both of these may depend on site; however, here we will 
assume that they are independent of site, though they may depend on lane. We say a 
model is lane-inhomogeneous (or simply inhomogeneous) if the transition rates and/or 
boundary conditions are not uniform between lanes; otherwise we say the model is 
(lane) homogeneous. A particle is blocked if there is another particle that prevents it 
from undergoing the forward motion below, otherwise it is unblocked. Transition rates 
may depend on whether a particle is blocked or not as illustrated in Figure [1] and the 
possible transitions we consider are defined below: 

• Motion: Plus-type particles move from (I, i) to (k,i + 1) with a change of lane / 
to k when unblocked (resp. blocked) at rate p 1 ^ (resp. p 1 ^). Similarly, minus- 
type particles move from {l,i) to (k,i — 1) at rate p l ^(b)- Motion is subject to an 

§ We identify lane M + 1 with lane 1 to represent this cylinder. 
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exclusion principle - a change can only occur if the target location to move into 
is vacant. Forward motion is used to mean the motion on a single lane (i.e., 
k — I) in either direction; this can occur only to a particle that is unblocked. If this 
forward motion is lane-homogeneous then we write 



• Particle type/direction change: Plus-type particles can change to minus-type 
particles with a change of lane / to lane k at rate p 1 ^. Minus-type particles can 
change to plus-type particles with a similar change in lane at rate p l Z^ (we assume 
that the site % is preserved for a change in type). For lane-homogeneous direction 
changes on the same lane, we write 



V- 



Boundary conditions for the model are assumed as follows: 

• Inflow: Plus-type (resp. minus-type) particles are injected at rate a l + (resp. a l _) 
into the minus end (reps, plus end) of the /— th lane. Both cases are subject to an 
exclusion principle. 

• Outflow: Plus-type (resp. minus-type) particles exit from the plus end (resp. 
minus end) of the /— th lane at rate f3 l + (resp. /3 l _). 
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Figure 1. Schematic diagram showing transition rates for the multi-lane bidirectional 
ASEP model; we consider M lanes on the surface of a cylinder. The plus end (resp. 
minus end) of the MT is at the right (resp. left). Plus-type (resp. minus-type) particles 
move step forward with rate p + (resp. p_) while they move forward associated with 
a change of lane with rates P^/^j (resp. p^S) when blocked (resp. unblocked). 
Plus-type (resp. minus-type) particles are injected into the left (resp. right) boundary 
of the system at rate a l + (resp. a l _) and exit at rate fi l + (resp. /3 l _) in the I— th lane. 
Particles can also change from plus-type to minus-type (resp. from minus-type to plus- 
type) with rate ptl (resp. p l Zt^) associated with a possible change of lane I to k (if 
I ^ fc); we assume the site is preserved during a change of type. 
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The state of the model {T l ±i {t)} at time t changes according to the above 
events (assumed to take place independently and instantaneously) from an initial state 
{7±j(0)}. In summary, the model for a given lattice (fixed M and N) has a number of 
parameters 

l-^k J->fc l^k l~^k I j nl 

P±,b 5 P±,u 1 P+- ) P-+ > a ± anQ Pi 

for k,l = 1, ■ ■ ■ M that collectively determine the bidirectional transport behaviour. We 
write a + = Z)|=i f a + to denote the total inflow of plus-type particles. In practice, 
many of these rates will be zero - for example, here we only permit a type change of 
a particle that remains on the same lane or moves to an adjacent lane. As in other 
ASEP models, we assume that the dynamics for typical choices of parameters converges 
to a unique statistical equilibrium independent of initial conditions, though possibly 
after a long transient. Exceptions occur if there is a degeneracy of the parameters, 
for example if no lane-changes are permitted, i.e., p l ±&u) = f° r k,l — ■ M, then 
accumulations can be created at arbitrary locations that are never cleared. The system 
exhibits various different dynamical regimes depending on parameters. We have not 
attempted to characterize the phase diagram of these states in full, but do discuss this 
further in Section fl~Tl 

The equilibrium densities of plus- and minus-type particles are defined to be the 
mean occupancy of the particles: 

d:=«i(*)>, *{:=««(*)> (1) 
where (•) denotes the ensemble average; assuming ergodicity this can be evaluated using 
a time average. In particular, if there is an accumulation at the tip then we denote the 
lane tip length \ l tip (t) as the number of particles in the accumulation within lane I at 
time t, and the tip size n tip (t) is 

M 

n t vp{t) = E A ^W- 
i=i 

with mean tip size as {n tip ) at a steady state. 

There are few methods giving exact analytical solutions of the equilibrium state 
for ASEP models (we refer to [27] for a review of ASEP models and their analytical 
solutions). Therefore, we numerically simulate this discrete event, continuous time 
model using both Gillespie |28j and fixed time-step (h t ) Monte Carlo methods. The latter 
method converges to the Gillespie method and gives outcomes that are independent of 
update methods in the limit ht — >■ 0. This general multi-lane model can reduce to 
previously studied ASEP models in special cases. For instance, restricting to M — 1 
and one particle-type gives the simplest unidirectional ASEP model. In particular of 
bidirectional transport, this multi-lane model limits the two models recalled in the 
following sections. 
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Table 1. The transition rates (a) with boundary conditions I in (b) for a two- lane 
model (M — 2) that gives bidirectional transport where each lane supports transport 
in one direction only and allows lane changes only together with a type change; if also 
associated with boundary Conditions II, the model corresponds to no-flux boundary 
at the plus end. 

(a) Transport rates 

p'_t=P^ = 0, p£S(g = 0(/ = l,2), p d :=pl-l 2 , Pu-=P^\ p± :=p£l 

(b) Boundary conditions I Boundary conditions II 
a\ = al = P\ = a 2 _ = 



2.1. A simple two-lane model 

A lattice model for bidirectional transport needs to either segregate, or deal with 
collisions between, opposite-directed particles. Works such as [12] include the 
possibility of binding/unbinding from/to a reservoir in their model. In our multi-lane 
model, for M = 2 one can segregate the particles into two lanes of unidirectional motion 
by setting the parameters in Tabled] (a) with boundary conditions I as in Tabled] (b). 
This particular choice of parameters corresponds to the two-lane model of [HI [19] and 
segregates the particles into two lanes- one lane will have no particles of minus-type, 
while the other lane will have no particles of plus- type, i.e., the equilibrium densities 
p 2 = cr 1 = 0. If, in addition, boundary conditions II in Tabled] (b) are satisfied then there 
is no flux at the right boundary and the model corresponds to the special case studied 
in [20]. We recall from [20] that in these circumstances, if the minus-end boundary 
condition a\ = a + small enough on the first lane and /3_ = P- on the second lane then 
a shock forms near the tip that traps a number of particles there. Using a mean-field 
approximation and denoting x = i5 where 5 = 1/N, one can find an asymptotic profile 
for the equilibrium densities of p and a for plus- and minus-type particles along the 
domain as follows (assuming a 2 and minjp 1 , (1 — p 1 )} are of order 5): 



a+ 



X 



(Npu _ Np d \ 



l + n£L(x-l), xe[x s ,l] 



x G [0, x s 



p+ 



[ ] ~ I i- P \x), xe[x s ,l] ' {6} 



where x s represents the shock position. Continuity of the flux J2U] then implies that the 
mean tip size can be approximated as 

2.2. A homogeneous thirteen-lane model 

Taking M = 13 and considering lane- homogeneous transition rates as in Table [2] (a,b), 
the multi-lane model reduces to the thirteen-lane model of [21]. If we consider boundary 
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conditions 

a l + = a+/M, a l _ = 0, l + = and l _ = p + (5) 

for all / G {1, • • • , M}, this corresponds to homogeneous injection of plus-type particles 
into the minus-end and no flux at the plus-end of the domain, while minus-type particles 
exit without impediment. Experiments and other information reported in [21J suggest 
the rates listed in Table [2] (a) are appropriate for this system. Rates in Table [2] and 
boundary conditions in (jSJ) are considered as default parameters for simulations unless 
otherwise specified. 

Table 2. Transition rates for the thirteen-lane model based on in vivo experiment 
measurements or estimations of velocities, mean run length and fluxes on a MT for 
the Ustilago maydis hyphal tip as detailed in [21] , Note that a+ represents the total 
injection rate, i.e., a + = YliLi a +- Simulations use N = 1250 sites (i.e., L = 10 fxm in 
length discretised with spatial step h s = 8 nm) and a sequential Monte-Carlo update 
with time step h t = 0.0042 s unless otherwise specified. 



(a) Possible transition rates s 1 


P+ P- P+- P-+ 
212.5 203.83 0.0406 0.0273 


P+.u 




J-M±l i-W±I 
P+.b P-.u 

4.335 


P-,b Q-+ f3 l _ 
106.25 1.06 212.5 


(b) Only adjacent lane changes 


permitted 


Type changes 


occur on a same lane 


p^)=o(MM±i) 




P-+ — P+- - 


= (fc^Z) 



3. Influence of collisions and lane-changes on transport properties 

In this section, we consider the multi-lane model with lane-homogeneous transition 
rates satisfying Table |2] (b) and boundary conditions as in (j5J) except for possibly 
inhomogeneous injection rates. For lane-homogeneous and symmetric lane-change rates, 
we remove the superscript and write 

P+,u(b) ■— P+, u (b)i P-,uQ>) --P-,u(b) 

independent of lane I. 

3.1. Mean field approximation and cross-lane diffusion 

In cases of inhomogeneous boundary conditions, the densities of plus- or minus-type 
particles can be homogenised by lane changes of particles, and lane changes are necessary 
for this to happen. Depending on which particles change lane under which circumstances 
(blocked or unblocked), this homogenization can occur to one or both types of particle. 

For example, if only minus-type particles can change lanes when unblocked in a 
dilute situation, then the density of minus-type becomes lane-homogeneous by cross-lane 
diffusion (see Figure [2] (c)), while the density of plus- type particles (not shown) may 
still be inhomogeneous due to the inhomogeneous injection rates and no lane changes of 
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Figure 2. Densities of plus- or minus- type particles corresponding to default 
parameters in Table [5] and boundary conditions in (|S|) except that the lane-change 
rates when unblocked are varied and there is injection only in one lane (namely 
a 7 + = 1 s^ 1 and a l + — for I ^ 7). (al) and (a2) show the inhomogeneous densities 
of plus- and minus- type particles respectively, when p± u = 0, i.e., particles only 
change lanes after collisions; (b) (resp. (c)) shows cross-lane diffusion leading to a 
more homogeneous density of plus-type (resp. minus-type) particles away from the 
boundaries with lane-change rates when unblocked as p~, u = 0, p+. u = 4.335 s _1 
(resp. p- >u = 0.4335 s , p+, u — 0). In all cases the densities near the tip (not shown) 
are high. 



this type. Even if only one type of particles changes lanes when unblocked, the density 
of both types of particles will be smoothed as shown in Figure [2] (b) for the density 
of particle-type particles. Exceptions to this are shown in Figure [2] (al,a2) when lane 
change only occurs after collision - in this case the plus-type particles injected into the 
middle lane "sweep" that lane clear of minus-type particles. 

These effects can be understood in terms of cross-lane diffusion of particles induced 



by lane changes in a dilute region. If we apply the mean field description of | Appendix A 



to a dilute region and assume small variation of density with site x = i8 {x G [0, 1] and 



5 = l/N), in the case of p_ >u = p +)U = 0, the solutions ( lA.lj) and (1A.2j) in Appendix A.l 



indicate the non-homogeneity of densities in both types of particles if given non- 
homogeneous injection rates a l + . If we allow lane changes when unblocked and assume 
that a change of type occurs more slowly than a change of lane, i.e., p~+,p+- P±,u, 
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then the stationary state distribution p l (x) satisfies (see Appendix A.l for details): 



n r dp 1 

= -5p + — +p H 

ax 



p 1 - 1 + p l+1 



2p l -5 



dp 



i-i 



dx 



8 



dp 



i+i 



dx 



1 = 1,2, 



M 



to leading order in 5. Solving this differential- difference equation for M lanes gives 

exp [— n + (k)x] 



LM/2J 
k=0 



J i 



2itkl „ 2nkl 
Ak cos — — — h B k sin 



where 



K + (k) 



M 



2 P+,u 



AI 



1-cos^ 



5p + + 2p + ,J cos 2nk ' 



M 
J i 



and Ak, Bk are constants determined by the density at p (xq). Similarly, the density for 
the minus-type particles can be approximated as 



a\x) 



LM/2J 

£ 

k=0 



2-wkl „, . 2nkl 
A k cos — 7^ — h iij, sin 



where 



K_(Jfe) 



M 



2p_ 



1 — cos 



M 



2nk 
M 



exp [k_(/c);c] 



5p_ + 2p_ iU 5cos^' 



and Ajj., B' k are constants determined by c 1 (xq). This suggests, not surprisingly, that the 
main effect of the lane-change in dilute cases is simply a diffusion of the density across 
the lanes in the direction of travel. 

Assuming there is a unique equilibrium state, the density for plus- and minus-type 
particles in the multi-lane model with homogeneous rates will have lane homogeneity 
though it will depend on site. For small injection rate a + , the system tends to a 
stationary state that is dilute near the minus end and in a high density near the plus 



end. Using a similar method as detailed in |Appendix A| to the lane-homogeneous case 
[20J, for small injection rate, under the conditions that p,a = 0(5) far away from the 
tip and ignoring 0(5 2 ), we find similar approximate expressions as (j2J) for the stationary 
state in dilute regions; see (1A.3j) . For large enough a + , dynamically pulsing states or 
filled states as discussed in Section 14.11 will appear. 



3.2. A phase transition in the tip accumulation 

For the half-closed system (i.e. closed at the plus end) in both the two-lane [20] and 
thirteen-lane models [21], it is clear that for typical particle- type change rates (e.g., 
P+- ^ P+)i particles accumulate at the tip. On the other hand, from the mean tip 
size approximation for the two-lane model, we see that (n t i P ) increase linearly with 
injection rate a + . However, this is not necessarily the case for the multi-lane (say, 
M = 13) model. In the following, we investigate the size of the tip accumulation and 
its dependence on the total injection rate a + for the thirteen-lane model with otherwise 
default parameters as in Table [2] and boundary conditions as in (jHJ). 
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Figure 3. (a) shows mean tip size as a function of total injection rate a+ for the 
multi-lane (M=13) model. Simulations PI refers to lane-change protocol PI where 
= p+/2 and p+.b = while simulations PII is for protocol PII where p±,b = P+/2- 
In both protocols we find a good fit to a rational function of a+ with a singularity 
at a critical value ot c ,i ~ 3.406 s _1 (fitting PI) or a C) jj rj 1.947 s _1 (fitting PII). For 
comparison, the straight dash line shows the linear dependence ((JJ) of a+ in the two- 
lane model recalled in Section 12. 1( other parameters use the values of corresponding 
rates in Table [2j see [20] for details, (b) shows the total occupancy J2i % T ±.i as a 
function of time for two typical injection rates as indicated, beyond and below the 
critical value for each of the lane-change protocols PI and PII. Other parameters used 
in the multi-lane model are default as in Table [2] and boundary conditions as in (|5[) . 



We compare two lane-change protocols PI and PII (further protocols are considered 
in the next section). The protocol PI only allows minus-type particles to change lanes 
when blocked while protocol PII allows both type to change lanes when blocked; lane- 
change rates are assumed homogeneous and symmetric (see details in Section 13. 31) . 
Simulations gave a mean tip size that depends on a + , fitting to a rational function 
with three fitting parameters A, B and C 

a + A(l + a + C) 



(n ti p(a + )) 



1 - a+B 



see Figure [3] (a) and Table [3] Note that there is a singularity in this rational 
fitting function at a + = 1/B. For the simulated system, these critical values are 
a C) i = 3.406 s _1 and a C) u = 1.947 s _1 respectively for lane-change protocol PI and 
protocol PII. Moreover, by comparing the change of the total occupancy of particles 
in the system in time between injection rates below and beyond the critical value, 
Figure [3] (b) indicates the existence of a transition between a shock state (where there 
is an approximately steady accumulation at the tip) and a new phase of motion, an 
unsteady pulsing state which is further examined in Section 14.11 

Additionally, by comparison with approximation (HJ) of the mean tip size from the 
two-lane model, we note that the mean tip size in the multi-lane model increases initially 
linearly but then nonlinearly for both lane-change protocols. This nonlinear increase 



Bidirectional transport and pulsing states in a multi-lane ASEP model 



11 



Table 3. Best fit of the data from simulations in Figure [3] (a) to a rational function 
(ntip(a+)) = a+ x~^~g ~ with parameters A, B, C. The values indicate best fit values 
and standard errors. 



A B C 



PI 
PII 



24.45± 0.329 0.2936± 0.00139 -0.2525± 0.00388 
22.18± 0.509 0.5136 ± 0.00908 -0.3978± 0.002219 



is due to trapped minus-type particles at the tip and is investigated in the following 
section. 

3.3. Influence of lane- changes on the tip accumulation 

The size of accumulation at the tip clearly depends on the lane-change rules; see 
Figure [3] (a). Recall from Section [37T1 that lane changes can homogenize the densities so 
we expect the homogeneity of X l ti (t) among the lanes will depend on the lane-change 
rules when blocked. To confirm this, we examine the maximum difference among the 
X l tip (t) defined as 



The distribution of this quantity gives a measure that characterizes how "smooth" (AX tip 
small) or "ragged" (AX Up large) the tip is. In particular, the larger the mean value 
(AXup), the easier it is for minus-type particles to be released from the tip. 

We consider four different lane-change (when blocked) protocols with homogeneous 
lane-change rates in an attempt to better understand their influence on the structure of 
tip accumulation via the distribution of AX(t) at stationary state. These are: 

(PI) = p+/2, p +) b = - only minus- type particles are allowed to change lanes with 
a rate (homogeneous and symmetric) that preserves the velocity; this corresponds 
to the assumption used in pi] ; 

(PII) p± ! b = p + /2 - both minus- and plus- type particles are allowed to change lanes; 

(PHI) p_6 = 0, p +t b = p+/2 - only plus-type particles are allowed to change lanes; 

(PIV) p l ^l~ x — P+/2, pJT& +1 = = p +: b - only minus- type particles are allowed to change 
and to only one of the adjacent lanes with lane- homogeneous rate. 

We ignore lane changes in the unblocked case as these are clearly not significant for 
escape from the tip. The distributions of the length differences AXu p at a fixed stationary 
time are shown in Figure H] for above four different lane-change protocols; other rates 
are as default ones. Note that allowing plus-type particles to undertake lane changes 
"smooths" the structure of the tip. The results for lane-change protocols PI and PII 
agree with the homogenized density of both type of particles on the effect of lane changes 
(see Figure [2]). Note the similar distribution of Au p between PI and PIV (PII and PHI) 
from Figure HJ it is then reasonable to consider protocols PI and PII as typical lane- 
change protocols. 



AX tip (t) = max{AL (t) - X k tip (t) : k, I = 1, • • • 



M}. 
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For both lane-change protocols PI and PII, particles can be trapped under several 
layers in the tip accumulation. Shown in Figure [3] (a) that the mean tip size increases 
more rapidly with injection rates for protocol PII, we suggest this is due to there being 
more minus-type particles "trapped" in the tip which in turn lead to larger increase 
on the tip size. We say a particle is "trapped" if its motion (including forward motion 
and lane changes) is obstructed by occupation of other particles and loss the ability to 
take the motion due to the exclusion principle. We define two fractions: the trapped 
minus-type fraction F_ )4 . rap and the minus-type fraction F_ within the tip in the steady 
state as 

mean number of "trapped" minus-type particles at the tip 
' trap mean number of total minus-type particles at the tip ' 

mean number of total minus-type particles at the tip 

( n tip) 

In the two-lane model discussed in Section 12.11 where two types of particles move on 
different lanes, the density profile of minus-type particles is dilute, which indicates 
F-^trap ~ 0. Additionally, from expressions (j2J) and (jSJ), 

a(x)dx/ / (p(x) + a(x))dx « - — (n Up ) 

JXs Zp + 

For the particular parameters in Table [2] with p^ = p H , the minus-type fraction 

F_ w 0.2% associated with a mean tip size (n t i P ) ~ 20. However, in the thirteen 
lane model, from Table HJ for low a + = 0.8 s^ 1 in lane-change protocol PI, F_ 4% 
associated with (n tip ) « 20 and F_ >trap ps 40%. This shows a difference between the 
two-lane model and the thirteen-lane model even for low injection rates, though with a 
similar mean tip size. Moreover, in the two-lane model, viewed as a queueing process, 
the expected delay for particles at the tip is while in the multi-lane model, the 

positive fraction of trapped minus-type indicates that the average time spent at the tip 
is larger then This average tip delay will be discussed later on in Section E33 

Comparing between the two lane-change protocols from Table HI protocol PI gives 
lower F_ tTap and F_ than those in protocol PII, with the same a + or the same mean tip 
size. This agrees with that a smoother structure of the tip accumulation is more likely 
to trap minus-type particles at the tip and indicates that the tip size and the (trapped) 
minus-type fraction are affected by each other. 



3.4- Measures of bidirectional transport efficiency 

The previous subsections indicate that allowing only minus-type particles to change 
lanes as a result of collisions is more effective (with less trapped minus-type particles) 
than allowing both types of particles to change lanes. To understand more on the 
bidirectional transport observed in living cells, it is interesting to speculate on the 
implications of the in vivo experiment data, particularly in Table [2j in terms of some 
non-trivial quantifiable efficiency 
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Figure 4. Distribution of maximum length difference AXa p for lane-change protocols 
PI-IV; other parameters are default as in Tableland boundary conditions as in ([5]). 
Note that the structure of the tip accumulation in PI or PII is clearly smoother than 
that in PI and PIV. 



Table 4. Comparison of trapped minus- type fraction F-^rap and minus- type fraction 
F_ at the tip between lane-change protocols PI and PII. This shows for protocol PI 



that there are lower F_ and F- 



.trap 



even for a similar mean tip size. When a + > oe c 



in both protocols, F-^rap ~ 100% and mean tip size does not convergent with time; 
marked as "— ". Other parameters are default as in Tableland boundary conditions 
as in ([5]). 



a+(s- 



0.8 



1 



1.2 



1.4 



1.6 



1. 



2.0 



PI 




37.8% 49.0% 
3.7% 4.5% 
19.7 24.9 



66.7% 

5.7% 
31.0 



72.5% 
7.6% 
36.8 



79.6% 
9.4% 
43.7 



85.7% 
12.1% 
52.1 



88.2% 99.7% 
14.3% 51.3% 
59.9 



PII 



59.7% 83.3% 
5.4% 9.9% 
20.3 28.1 



89.8% 
14.3% 
36.5 



92.6% 
20.2% 
48.3 



95.8% 
28.6% 
70.1 



98.5% 
42.7% 



99.6% 100% 
50.4% 57.6% 



There are many possible ways to quantify efficiency of bidirectional transport, 
depending on what sort of behaviour is important. For example, in some circumstances 
it could be that maximum speed is required while in other circumstances it could be 
that the maximum flux is required - for a practical system of bidirectional transport 
on a MT, different notions may be needed for different purposes. In a fungal model 
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system, the particles accumulate at the MT plus end are thought to prevent endosome 
falling off the MT and to support their backward transport [21]. Therefore, efficiency 
transport of dynein requires a certain proportion of dynein reaching the tip to form an 
accumulation and a short duration at the tip. Therefore, we consider particle arrival 
efficiency and average tip delay defined as following. 

Suppose that the n-th particle injected at the minus end is at location = 
(L n (t) , I n (t)) at time t, and suppose that: 

• it enters the system at time t n ,\ and leaves the system at time t nA ; 

• the farthest point it reaches in the system is d n , i.e. d n = maxt /„(£); 

• it enters the "tip" at time t n>2 , i.e. this is the smallest t such that I n (t) > L tip ; 

• it leaves the "tip" at time t n ^, i.e. this is the largest t such that I n (t) > L tip ; 

where L tip refers to the farthest site away from the plus end at the tip, typically we set 
Lup ~ 90%iV for simulations. These times are illustrated in Figure [5] (left panel) with a 
comparison of endosome motility in vivo in Figure [5] (right panel). If the particle does 
not reach the tip then we have only t nt \ < t nA defined, while if d n > Lu p then we have 
t n ,i < t n ,2 < £n,3 < ^n,4 defined. Based on these notations, we define the particle arrival 
efficiency to be [Jj] 

E 1 = lim E 1 (T) = lim f r {w ' £ ~ ^ ~ < T \ ■ 

Note that < E\ < 1 is a dimensionless measurement of efficiency that is close to zero if 
few of particles arrive at the tip before leaving the system while it is close to one if most 
of particles arrive at the tip before leaving the system. The quantity E\ is independent 
of how long it takes to arrive at the tip or leave the system. Meanwhile, we define the 
average tip delay to be 

E 2 = lim Eo(T) = lim ^ : \\ dn \\>L Up ,t nA <T}(tn,3 ~ t n , 2 ) 

r->oo 2V ; t->oo #{n : d n > L tip , t nA < T} 

The average tip delay E 2 has unit of time, and measures the average time spent at the 
tip for all particles that get there: larger values indicate that particles are trapped at the 
tip for a long time while small values indicate that the particles wait only a short time 
before leaving the tip. Note that the quantity E 2 is independent of E\, the proportion of 
particles that arrive at the tip. Note also that if we assume homogeneous particle-type 
change rates then E 2 > l/p+_. 

One can derive approximate upper and lower bounds of E\ for the lane homogeneous 
multi-lane model. Considering the density of plus-type particles in dilute situation; see 
density approximation (1A.3I) in Appendix A. 2 , we have 

E x < min{p(l)/p(0),l} (6) 

Np_ + Np + _ 



mm < exp 



P- + P-,u i P+,u P+ + P+,u i P+,u 



1 



T is assumed to be large enough that any transients have decayed, i.e., t\ i > T . 
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Up Plus end plus end 




Figure 5. Typical motion of three particles injected into the system (left panel). 
We show the times tj :in m = 1, ...,4 for the j—th particle. Particles that exceed 
L t ip are said to have entered the tip region. Kymograph (right panel) shows motility 
of early endosomes that were visualized by photo-activation of the endosome marker 
paGFP-Rab5a [21] in subapical regions of a cell of the fungus Ustilago maydis. Note 
the similarities in the two panels on the trajectories of motions. 



while considering the probability of first particle-type change occurring at the tip, we 
have 



Ei> 



P+ + P+,u + P+,u 



1 ^ (7) 

as particles may change directions (type changes) several times before reaching the tip. 
Figure El (a) shows the upper and lower bounds for E\ compared to the simulations. 
Moreover, from Figure El (a) and Table this arrival efficiency is independent on the 
particle-type change rate. 

The average tip delay E 2 is clearly dependent on p + _ for a homogeneous model as 
l/p+- is the expected delay of minus-type particles without obstruction such as in the 
two-lane model recalled in Section [2TT1 (see [20] for details). We therefore consider the 
delay ratio E 2 X p+-, i.e., the ratio of the average tip delay to the expected delay if there 
are no trapped minus- type particles. This ratio is usually no less than 1, and when it 
is close to 1 this is the most efficient situation in terms of the delay being minimal. 
Simulation results in Table El and the inset to Figure El show that for both lane-change 
protocols in the multi-lane model, this delay ratio is larger than 1 or E 2 > This 
agrees with results in Table H] showing the existence of certain proportion of minus- 
type particles being "trapped" at the tip. Comparing Figure El (a) and the inset in 
Figure El (b), we consider the dimensionless quantity E 2 x p + ^/E± to be a balance 
between arrival efficiency and delay ratio. As varying p + _ we find a minimum value 
where p + _ is closed to the experiment data of p + _ = 0.0406 s -1 ; see Figure [6] (b). This 
offers an implication of in vivo transport rates. 

Comparing the two lane-change protocols PI and PII, protocol PI is more effective 
in the sense that it has a relatively lower delay ratio E 2 x p + ^ when varying a + (see 
Table [5]) or for small p + _ (see the inset in Figure El (b)). Similar to the case for trapped 
minus-type fraction and mean tip size, the average tip delay time increases more rapidly 
with injection rate in protocol PII than the increase in protocol PI. 
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Figure 6. Illustration of E\ and E2 when varying (a) shows that the arrival 

efficiency E\ in both lane-change protocols PI and PII behaves similarly as varying 

p^ with approximate upper bound from ([6]) and lower bound from ([7]). (b) shows 

the balance between delay ratio and arrival efficiency quantified as E2 x p^ / ' E\ for 

protocol PI; the inset shows the comparison between two protocols PI and PII of the 

delay ratio Ei x p^ ; errorbars are indicated. Simulations are done by averaging over 

10 runs using the Gillespie algorithm in a time interval [To,T] = [500s, 3000s]; other 
rates are default as in Table [2] and boundary conditions <j5j> . 



Table 5. Comparison of E\, E2 between lane-change protocols PI and PII. The arrival 
efficiency E\ in both protocols is relatively independent of a+ and drops to almost zero 
when a + > a c j for protocols PI (a+ > a c ji for protocol PII). The average tip delay 

E2 in both protocols are larger than the expected delay l/p_| « 24.6 s~ 1 when no 

minus-type particles are trapped, and the delay time increases more rapidly for protocol 

PII than for PI when increasing a+. Note that the particle- type change rate p^ is 

fixed here, meaning that the delay ratio E2 x p^ is a multiple of E^- Simulations are 

done as in Figure [B] with otherwise default parameters. 



a + (s- 1 ) 0.8 



0.9 



1 



1.1 



1.2 



1.3 



1.4 



PI 



Ei 

E 2 (s) 



81.9% 81.9% 
30.71 30.87 



81.5% 
31.1 



81.9% 
31.1 



81.6% 
31.8 



81.6% 
32.6 



81.9% 
33.5 



81.7% 



PII 



£1 

E2 (s) 



81.4% 
31.9 



81.8% 
33.2 



81.7% 
34.2 



81.7% 
37.9 



81.9% 
44.9 



81.8% 
48.7 



81.5% 
61.4 



4. Critical behaviour and pulsing states 

Results from Section 13.21 suggest that the mean tip size does not converge when an 
injection rate a + approaches a critical injection rate a c j for lane-change protocol PI 
(a c ji for protocol PII) in a finite system; additionally, by simulating the efficiency 

defined in terms of delay ratio E 2 X p^ , the system will be less efficient with either a 

large total injection rate a + (see Table [5]) or a small particle-type change rate p + _ (see 
Figure |6] (b)). Large mean tip size or inefficiency on E% x p H are expected due to high 
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densities of plus- and minus-type particles, and in the following we discuss the dynamics 
of pulsing states that exist in this model. 

4-1. Pulsing states in the system 

For injection rate a + smaller than the critical value a c discussed in Section I3.2[ the 
system converges to a stationary state (shock state) with an accumulation of certain size 
at the tip, and this tip size increases as a+ approaching a c . For a+ > a c and otherwise 
default parameters, a new type of behaviour that we call a pulsing state can appear. 
These are characterized by the behaviour in Figure [7J For finite systems, near a c the 
accumulation at the tip may show large irregular oscillations including a propagating 
region or pulse of high density that moves away from the tip before dispersing. 

In a pulsing state, the whole system density undertakes large approximately periodic 
oscillations. These can be split into two parts: 

• A filling phase during which the density at the minus end is low, and a growing 
high density pulse of mixed particles moves steadily towards the minus end. During 
this phase there is a net flux into the system. 

• An emptying phase during which the density at the minus-end is high, and the 
pulse propagates out of the system. During this phase there is a net flux out of the 
system. 

For most of the cycle there is a low density region between the tip and the pulse. 
The pulsing occurs through an alternation between these phases as fronts (evident in 
Figure [7J separating low and high density regions move through the domain. At other 
parameters we find the system converges to a filled state for a + larger than a c . Figure El 
shows examples of the pulsing and filled state. In the filled state, once the filling phase 
is complete the system remains approximately full. 

To better understand the behaviour of system in a pulsing state, Figure [9] shows 
successive horizontal lines as the evolution of the medium time average of the density 
Y^i=i{ T -,i + T +,i)/M (over successive blocks of 1260 s) from a time-series and gradual 
changes of the local density. Figure |9] (a,b) shows rapid convergence to a stable 
accumulation at the tip, though (b) displays more fluctuations. Examples of pulsing 
states are illustrated in Figure [9] (d,e,f ) for increasing values of a + ; these exhibit long 
approximately periodic oscillations. Figure |9] (c) shows an intermediate situation where 
the tip accumulation shows irregular large amplitude oscillations. The main behaviour 
in a high-density mixed region will be changes of particle types and therefore we expect 
the ratio of densities between minus- and plus-type particles to reach the equilibrium 
ratio given by the type-change rates, i.e. 

P-+ ~ P+- 
p = k, a = k 

P-+ + P+- P-+ + P+- 

where k = p + b « 1 represents the total local density. 
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Figure 7. Time progress of a pulsing state for the thirteen lane model of size N = 100, 
a + = 4.0 and no flux on the plus-end; otherwise parameters as in Tableland 
boundary conditions as in ([5]). The sequence of filling (top three frames) and emptying 
(bottom three frames) repeats approximately periodically. Observe that the pulse of 
high-density mixed particles propagates slowly to the left via diffusion of vacancies 
through the pulse. Grey indicates plus-type particles (that move to the right) while 
black indicates minus-type particles (that move to the left). The white regions indicate 
vacancies. 



For both lane-change protocols PI and PII we find pulsing states for a + > a c and 
otherwise default parameters - we interpret this as an existence of a negative net flux 
for the high-density mixed phase. Note that densities p and a are primarily governed 
by the particle- type change rates p + _ and if we increase then the net flux 
may change to be positive and the pulsing state will be replaced by approach to a 
high-density uniform filled state. For example, Figure [8] (a) shows the pulsing state for 
p_ + = 0.05 with a + = 3 while Figure El(b) shows a filled state for p_ + = 0.14 s _1 
and a + = 4 (other parameters are the default values). For the latter case, we infer 
that an initial pulse with j > will saturate to p + a = 1 and then j — > as the number 
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Figure 8. (a) shows a pulsing state for p |_ = 0.05 and a + = 3 s _1 ; (b) shows 

a filled state for p |_ = 0.14 s _1 and a + =4 and this approaches a homogeneous 

high-density mixed equilibrium state. The other rates are default as in Table [2] and 
boundary conditions in ([5]). 



of vacancies in the high density state goes to zero. It will be an interesting challenge to 
understand the statistical properties of these high-density mixed states and, for example, 
to predict j from the parameters of the system, as there might be more than one particle 
compete for a common empty site. 

Numerical simulations show that this pulsing state is robust to changes in 
parameters - Figure [HJ] shows phase diagrams for lane-change protocols PI and PII 
on varying a + and p + _ with other default parameters. The particle-type change rates 
we usually use obeys < which means that a typical switching of plus- type 

particle occurs at the tip (see E\ in Table EJ), this novel pulsing state also appears when 
E+^^E=± > i. see Figure [TTJ A necessary condition for a pulsing state we suggest is 
that a + is beyond a critical value (i.e., not in a shock state) and there is a net negative 
flux in the high-density mixed region. 

4-2. System size effects 

Although the previous section focuses on a fixed system size N, the critical transition to 
pulsing states is invariant of N as iV — > oo under certain scalings of parameters for the 
multi-lane model with boundary conditions (jSJ). For example, if we scale the forward 
motion rates, keeping p±/N constant and use the lane change rates when blocked p+/2 
as default whilst keeping the particle-type change rates constant, we find no significant 
variation of mean tip size with iV as illustrated in Table El This agrees with the 
approximation from the two-lane model. For injection rates beyond a critical value 
and the same scaling we find pulsing states for arbitrarily large N. 
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Figure 9. Each horizontal line represent a time-average of the occupancy Yli=i ( T - i + 
t 1 + j)/M over blocks of length 1260 s in the homogeneous thirteen-lane model (M = 13) 
with injection rate a + indicated (other rates as default). In (a), the model with a 
small total injection rate reaches its stationary state with a small number of particles 
accumulated in the tip. (b)-(d) show the density-time courses for injection rates near 
the critical value. In (d)-(f), pulsing states appear where the accumulation moves away 
from the plus end and simultaneously grows. 



Table 6. The tip size when varying the system size N for parameters in Table [2] 
and boundary conditions ([5]) except keeping p+/N — 0.17, p~/N = 0.1631 and 
p_.b = p+/2. There values are chosen to be consistent with rates in Table [2] when 
N = 1250. Mean and standard error of mean are shown. 



N 
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Figure 10. Numerically obtained phase diagrams for lane change protocol PI (left 
panel) and lane change protocol PII (right panel). Other default parameters are as in 
Table [2] and boundary conditions as in ([5]). Bars indicate uncertainty of the borders 
between different states. 




4-3. A simple model with a critical injection rate 

To better understand the appearance of a critical injection rate a c that leads to a 
transition in the behaviour of the tip, we introduce a simpler model that computably 
predicts nonlinear behaviour of the mean tip size with injection rate a+, and a singularity 
of the mean tip size at finite injection rate. 

Consider a single lane first- in, last-out queue of plus-type and minus-type particles. 
Plus-type particles are assumed to arrive at the left end of the queue at a rate a per 
time-step. Every time-step, we assume that all particles independently to be minus- 
type with probability p and plus-type with probability 1 — p. A number of particles is 
assumed to leave the queue whenever the left-most particle is of minus- type, in which 
case all of the adjacent minus-type particles are assumed (instantaneously) to leave the 
queue. This results in a steady growth and intermittent loss of particles. At a given 
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time-step t, if we assume the queue has size n(t) then after random changes in type, the 
number of particles that leave the queue can be expressed as 

n \—p n 

^( n ) = m p m (i — p) + n p n — p~ 



m=l 



1 — p 



Hence the mean rate of net growth of the queue will be given by a + — A(n). This will 
give zero net growth when a + — A(n) = 0. Therefore, with injection rate a> + , the mean 
size of the queue is 

Infl-aM 

("(«+)> = -h 

lap 

which predicts a critical injection rate a c = For a + < a c , (n(a + )) grows 

monotonically and nonlinearly with lim a+ ^ Q , c _ (n) (a) = oo. These are clear analogies 
with the observations of mean tip size in the thirteen- lane model in Figure |3] (a). Indeed, 
one can fit the data in Figure |3] (a) to a logarithmic function but the fit is inferior to the 
rational function discussed there, suggesting that the model above is too simple for an 
accurate quantitative explanation. This simple model reproduces the nonlinear increase 
for a + < a c in Figure [3] (a) and the singular behaviour at the critical value a c . For 
a + > a c it predicts unbounded growth of the queue, analogous to convergence of the 
system to a "filled state". 



5. Discussion 



Motor-driven bidirectional transport of vesicles and organelles is vital for the 
organization and function of eukaryotic cells and intracellular motility serves various 
cellular processes [29]. The basic function of bidirectional vesicle transport is to deliver 
cargo over distances, thereby modifying gradients and ensuring communication between 
different regions of the cell. 

By considering the process in a number of stages, starting with cargo uptake, 
followed by transport along the fibres of the cytoskeleton, and finally ending with cargo 
off-loading, the site of unloading is often also a region of cargo uptake. Transport back 
on the cytoskeletal track therefore not only recycles the transport vesicles, but also 
serves for long-distance delivery to other regions of the cell. In a fungal model system, 
it has recently been shown that dynein are concentrated at MT plus-ends to prevent 
organelles falling off the track and this concentration of dynein is done by an active 
retention mechanism (based on controlled protein-protein interaction) and stochastic 
motility of motors [21, 26J. 

An accumulation of motors at MT plus-ends can be seen as an inefficiency, as 
dynein motors are supposed to transport particles rather than wait at the MT end. 
Thus, the cell needs to find a compromise between (a) ensuring that organelles are 
captured at MT ends, which requires dynein accumulation at the plus-ends [2T] and 
(b) keeping dynein moving along MTs to deliver the cargo to minus-ends. As discussed 
in Section 13.41 an interesting result of this study is that the parameters from [2T] in a 



Bidirectional transport and pulsing states in a multi-lane ASEP model 



23 



Ustilago maydis system (and p + _ in particular) do address this compromise between (a) 
and (b). Moreover on comparing different lane-change protocols (that is, either dynein 
changes lanes or both kinesin-1 and dynein change lanes), we find that it is necessary 
for only dynein to change lanes in order to give biologically realistic system states over 
a larger range of possible fluxes; see Figure [3] (a). 

The maximum unidirectional flux in the thirteen-lane model can occur, by 
analogy with the so-called maximal current state in the ASEP model [3], [27], when 
a + ~ 13 x p+/2 « 1300 assuming appropriate boundary conditions and only 
one particle type. Note that the maximal flux with bidirectional transport in the 
thirteen-lane model is about the half of that in unidirectional transport. However, 
the existence of a critical injection rate a c <C 1300/2 s _1 (see Figure [3]) in the half- 
closed homogeneous thirteen-lane model places a much lower bound on the maximum 
flux. The experimentally measured a + ~ 1.06 in Table [2] is of the same order as 
a c . To achieve bidirectional flux that is half of the maximum unidirectional flux (about 
600 s" 1 ) on a MT, much more organized transport with specific control of lane-change 
to segregate different particles into different lanes (as in Table [1]) is therefore necessary. 

For the half-closed system, we have found for a + > a c some novel "pulsing states" 
of the system. It will be very interesting to further explore the region of existence of 
these states. There may exist other new phases for other boundary conditions and/or 
transition rates in the multi-lane model. An open question includes whether one can 
find a better understanding of the dynamics of pulsing states and/or whether there are 
analytic approximations that confirm the phase diagram approximated numerically in 
Figure [TOl 

The model can be generalized to incorporate a number of features that may be 
important for MT transport. This includes (a) additional species of motor with different 
transport properties, such as the several types of kinesin known to exist in vivo pfl 30J; 
additional motors can be expected to increase collisions along the MT length and cause 
high density of motors that might affect the transport efficiencies described in Secton l3T4l 
(b) The detachment and reattachment of motors from the MT [TT] may be important 
in other systems; the size of any accumulation will clearly be affected by this, (c) 
Non-trivial geometry of MT bundles needs consideration as motors may jump between 
different MTs; this might be another possibility to avoid collision between counter- 
moving organelle, (d) A more realistic motor and cargo size will influence the ability 
of lane changes to overcome blockages. Finally (e) static obstructions such as MAPs 
(microtubule associated proteins) will influence the behaviour of motors on the MT by 
increasing more potential blockages j24j [3TJ [32J [33] ; a high concentration of these will 
clearly influence the transport efficiency At present, quantitative data for most of these 
processes are not available for any living cell. With each additional process one can 
gain quantitatively more precise mathematical representation of the full MT-mediated 
transport within a cell. It remains to be seen whether new qualitative effects and better 
understanding will result from more accurate models. 
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Appendix A. Mean field approximation for the multi-lane model 

We use a mean field approximation to give an equation for evolution of the densities 
p\ and a\ for the multi-lane model with lane- homogeneous rates (but not necessarily 
boundary conditions). Examining the balance of incoming and outgoing particles to a 
given site, the mean field approximation gives the following, where on the right hand 
side we write pi to mean p\, etc. 

dp\ \ /-< \ 

= P+Pi-li 1 - Pi- <Ti) +P- + 0~i -P+Pi{l ~ Pi+l ~ (Ti+i) ~P+-Pi 

+ plj-\{l - p\- 1 - <J^ l ){\ - Pi - a t ) 
+ pt,uP?-\(l ~ P? 1 ~ ~ Pi ~ *i) 

+ {plMOt 1 + °n + pi,ptM +1 + a - ^ - 

- (pIA 1 - p\~+\ - + pi A 1 - p l i+\ - Pii 1 - pr+i - Vi+l) 



There is a similar expression for but for reasons of space we do not give this here. 
As the general mean field equations are not easy to solve, we consider below two special 
cases of these mean field equations. 



Appendix A. 1. Dilute lane-inhomogeneous densities 
Let us assume that 

• p\{t) = p l (x,t) with x — iS (where 5 = 1/N = h s /L is small) and the spatial 
dimension is parametrized by x G [0,1], 

• p and a are small (i.e. dilute limit); 

• the lane changes rates are lane homogeneous and symmetric (i.e., p± >u = p±Ju); 

On expanding p\ +1 = p\x) + 5^ + 0(5 2 ) and discarding any terms that are quadratic 
in 5, the mean field equations above simplify to 

¥ r dp 1 ( l+1 o , dp 1 - 1 dp l+l \ 

+ p-W -p+-p 1 

da 1 Ba l ( l+1 n , da 1 ' 1 ,da l+1 \ 

— = 5p_—+ a l+1 + a 1 ' 1 - 2a l + 5—- + 5—- 
at ox \ ox ox J 

+ p+-p l -P-+0- 1 

which can be used to characterize the combination of bidirectional transport (p±), 

change in direction (p+_ and p h ) and cross-lane diffusion (p±, u ) on the density. 

Considering a region of the domain where the dilute approximation holds, the stationary 
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state distribution will therefore satisfy 



ax \ ax ax 



+ P-+o -P+-P 



. da 1 ( i , , n , .da 1 1 da l+1 \ , , 

-6p-— = p-, u a l+1 + a 1 ' 1 - 2a l + 5—— + 5—— + p+J - p_ + o l . 
dx \ ax dx J 

If p± u = 0, adding the above two equations gives p+^r = P-^ri meaning that in 
dilute situation, p and a have a linear relationship. This gives a solution to the above 
ODEs: 

" 'p~+N p+-N^ 



p{x) = C exp 
p + C 



x 



P- 



cxp 



V- P+ ) 

(p-+N p + _N 



+ D 



+ 



P+-D 
P-+ 



(A.l) 
(A.2) 



P- P^ 

where C, D can be determined by the densities at p(xo) and a(xo). 

If p± tU are large (p + __,p_ + <C P±, u ), ignoring the turning rates, the above ODEs 
lead to 

dp 1 
dx 
da 1 
dx 



<5Ph 



J + i , J-i oJ x^" 1 A.^ +r 



p+, M p (+1 + p- 1 - 2p< - tf- 



-5p_— = p_ tt <r m + a'" 1 - 2a ; + 5 



dx 
da 1 ' 1 



dx 



+ 5 



dx J 
da l+v 



dx 



Appendix A.2. Lane-homogeneous densities 

Now consider homogeneous model which gives lane-homogeneous densities p(x) = p\ 
due to the symmetric structure of the cylinder and as before x = id. In this case we 
have, ignoring second and higher order terms in 5 = 1/N, that the stationary state in 
the dilute case satisfies 

dp 



0=p. 







_<7 



Ph 



-P+ - P+ 



-P + S 

p-+o + p+_ P + 5 \ P _ + p!^ +i + jHr 1 



dx 
dp 

dx 



Assuming zero net flux, this implies that v+p — v_a = where v±/h s = p± 
p±w -1 - Applying the boundary conditions p(0) 



gives 



p(x) 



O: i 



p + M 



exp 



iVP-H 



iVp+_ 



P+ + P+,u + P+,u , 



P±t +1 + 



(A.3) 



„P-+P-,u +P- 

analogous to the expression found for the two lane model in |20j but spread over all M 
lanes. This mean field approximation of density works well for low densities but near 
the tip, plus-type and minus-type particles do not have complementary density as in 
[20] ; in fact they can be at the same order near the tip for typical turning rates; for 
this reason it does not appear to be easy to obtain the mean tip size from a mean field 
approximation. 



